Search graph structure and its implications for multi-graph constrained routing and scheduling problems

Multi-graphs where several edges connect a pair of nodes are an important modelling approach for many real-world optimisation problems. The multi-graph structure is often based on infrastructure and available connections between nodes. In this study, we conduct case studies for a special type of constrained routing and scheduling problems. Using the airport ground movement problem as an example, we analyse how the number of parallel edges and their costs in multi-graph structure influence the quality of obtained solutions found by the routing algorithm. The results show that the number of parallel edges not only affects the computational complexity but also the number of trade-off solutions and the quality of the found solutions. An indicator is further proposed which can estimate when the multi-graph would benefit from a higher number of parallel edges. Furthermore, we show that including edges with dominated costs in the multi-graph can also improve the results in the presence of time window constraints. The findings pave the way to an informed approach to multi-graph creation for similar problems based on multi-graphs.

Many optimisation problems in transportation, logistics or telecommunications can be formulated as search on a multi-graph. An example of problems include the vehicle routing problem, hazardous material transportation, multimodal shortest path problem and airport ground movement problem to name a few. The multiple parallel edges between pairs of nodes of the multi-graph offer a convenient way of modelling the real-world structure and inherent multi-objective nature of the problems including time, economic or environmental objectives. The parallel edges can represent routes with different costs in the multi-objective vehicle routing problem 1-3 and hazardous material transportation 4 , different modes of transport 5,6 in the multimodal shortest path problem and tour planning 7 , or different speed profiles in the trajectory based traffic management [8][9][10] . Furthermore, there are often various constraints which have to be satisfied by the solutions in order to be feasible, for example a delivery vehicle must visit customers in specified time windows.
So far, the research on the abovementioned optimisation problems focused mainly on the search algorithms for finding the best solutions using the multi-graph formulation of the problem. The structure of the multi-graph has been considered given and fixed, representing the real-world connections between the nodes. This is to some extent true for some problems such as the vehicle routing or multimodal shortest path problem where the multi-graph represents the underlying and existing infrastructure (roads, rail lines, etc.). However, even for these problems, infrastructure and schedule design is worth of investigation. For other problems such as the trajectory based traffic management, the multi-graph structure is mainly a result of the modelling approach, where the number of the parallel edges and their costs are design parameters. For example in 10 , there is an infinite number of speed profiles between two nodes (corresponding to different continuous speeds) and the multi-graph can include only some of them. As noted in the previous research 11 , the number of the parallel edges and their costs not only affect the computational times of the search algorithms but also the quality of the obtained solutions.
In this study, we further analyse the cases and conditions when the multi-graph structure has a direct consequence for the search algorithm. The case study of the airport ground movement problem is used to demonstrate how the decision on the number of the parallel edges and their costs affect the quality of the found solutions. The results in turn can inform the creation of the multi-graph not only for the airport ground movement problem but also other problems such as the above mentioned vehicle routing problem, hazardous material transportation, the multimodal shortest path problem and trajectory based traffic management. For example, the need to determine the number of multi-objective speed profiles for the trajectory based traffic management [8][9][10] and the vehicle routing problem including energy-efficient driving 6,9 plays a critical role in providing solutions, as it does in the airport ground movement.
The airport ground movement problem is a combined routing and scheduling problem which aims to find conflict-free routes and schedules for all aircraft taxiing between gates/stands and runway or vice versa with minimum taxi time and fuel consumption. The airport taxiway layout corresponds to a simple graph, which is expanded into a multi-graph by considering segments (i.e. a sequence of edges of the same type such as straight and turning) and their associated speed profiles as shown in Fig. 1.
The routing algorithm uses the multi-graph for the search and has to select: (1) which segments to include in the route between the start and end nodes; (2) which speed profiles to use for the selected segments. In order to prevent conflicts, each edge can be occupied only when it or a nearby edge is not traversed by another aircraft. For this purpose, aicraft must satisfy a time constraint called time window on each edge of its route. In this paper, we consider two objectives: obj 1 is the taxi time and obj 2 is the fuel consumption. A specialised routing algorithm AMOA* 11 based on multi-objective A* algorithm is used in this paper. It should be noted, that any search algorithm, e.g. metaheuristic 12 , can be employed for this purpose.
The routing algorithm considers aircraft iteratively according to their start times. This first-come-first-served sequencing has an advantage to consider aircraft sequentially as they become ready to start taxiing. For each aircraft, a set of routes with nondominated costs are found by the algorithm using edges with feasible time windows. The edges of infeasible time windows are avoided by the routing algorithm, causing a detour. If no route is found due to the lack of available time windows, the start time of the aircraft is postponed by 60 s. This value is set approximately as airports usually operate (e.g. estimated time of departure) with a precision in minutes 11 . The start time is iteratively extended until time windows become eventually available. Once a set of routes is obtained, one route with the minimum cost is selected. The minimum cost is calculated by multiplying the values of obj 1 and obj 2 with corresponding unit costs w P = (0.469, 0.71) . The unit cost of 0.469 EUR/s for taxi time was calculated in 13 and includes the cost for maintenance, fleet and crew. The unit cost for fuel consumed is set to 0.71 EUR/kg as in 13 . It should be noted, that w P is not utilised within AMOA*. w P used here replaces a decision maker who would in real operation select a route for each aicraft from a set of Pareto routes according to his/her preferences. After the route is selected, time windows are updated for edges belonging to its segments. To ensure a safe separation from other aircraft, also time windows of edges within a threshold distance of 60 m are blocked in addition to the edges of the selected route. The separation of 60 m corresponds to approximately 12 s difference between successive aircraft at taxiing speed 10 knots, similarly as in 13 . Table 1 summarises the definitions used in this paper.

Results
Airport ground movement problem instances. In this paper, we use a set of real instances of arrival and departure flights from 3 airports: Doha International Airport (DOH), Hong Kong International Airport (HKG) and Beijing Capital International (PEK). The complexity of the taxiway layout ranges from simple (DOH), medium (HKG) to complex (PEK), as shown in Fig. 2. The graphs and flights are detailed in Tables 2 and 3 . The data specifies landing/pushback times, gates/runway exits and the weight category for each flight. The instances contain traffic data as follows. Instances marked original are from our previous work 11 . The instance ins1 is an artificial instance with introduced aircraft conflicts. ins2 is similar to ins1 but with conflicting aircraft isolated such that the aircraft later in the sequence are not affected by routes of previous aircraft. Finally, instances with increased traffic (0-100 %) are used. The instances with introduced aircraft conflicts and higher traffic are expected to have less available time windows which may benefit from a multi-graph structure with higher u. The process of creating the airport layout and capturing the flight data is detailed in the "Methods" section. Besides the number of nodes in the multi-graph (depending on the airport), the number of parallel edges corresponding to the speed profiles (the size of C n,m ) significantly affects the size of the search space. As mentioned above, there is a large number of speed profiles (in this study, we used 20) available for each segment. As a result, the complexity rises exponentially as shown in Fig. 3 and experiments with more than 5 speed profiles for each segment could not be finished in reasonable time for all airport instances. To address this problem, two different speed profile selection approaches are used to reduce the multigraph. Figure 4a,b illustrate the difference of speed profile selection for increasing number of u for a route of a single aircraft. The number of found Pareto optimal solutions gradually grows with the larger value of u. With Table 1. Notations used throughout the paper.
The directed graph of airport taxiways with nodes n ∈ V and edges e ∈ E The vector of preferences (i.e. weights for objectives) representing unit costs  www.nature.com/scientificreports/ evenly distributed speed profiles, the found solutions gradually cover the whole Pareto front with a larger u. The same is true for the preference-based selection of speed profiles, where the solutions concentrate at the preferred region with a smaller u and gradually spread to cover the whole Pareto front with a larger u. Both speed profile selection approaches are described in detail in the "Methods" section.    In order to analyse the effect of u on the quality of the obtained results, a baseline with u = 1 is established. The multi-graph is reduced to a single-graph with a single speed profile for each segment. The preference w P = (0.469, 0.71) is used to select the speed profile and for reserving a route for each aircraft. Further experiments were conducted with u = 3 and u = 5 , where the speed profiles for each segment were selected evenly. In order to facilitate the comparison, Table 4 shows savings of total cost (calculated as obj 1 · 0.469 + obj 2 · 0.71 ) obtained with u = 3 and u = 5 compared to the baseline case. The values of objectives are in "Appendix" section. Positive values in column cost refer to the cost saving. The column marked as 'improved' refers to the number of aircraft better in both objectives with respect to u = 1 . In order to estimate the number of aircraft that can be improved, indicators I 1 and I 0 are proposed with more details in "Methods" section. I 1 is the number of routes which can be potentially improved using alternative speed profiles. I 0 is the number of routes which cannot be improved even by using alternative speed profiles. Negative values of cost refer to higher costs than those of the baseline case with u = 1 . As more speed profiles are included in the multi-graph with u > 1 than the baseline case, finding worse routes should not be possible. Worse routes are caused by sequential routing of aircraft. In individual cases, the routing algorithm searching the multi-graph with larger u can always find a route with better costs. However, this can have detrimental effect on the availability of edges for the subsequent aircraft, resulting in a detour and higher objective values.
To gain more insight, firstly, we analyse the instances marked original. The higher values of u resulted in mostly negative savings < 1% which are caused by sequential routing of aircraft as explained above. The lack of improvement can be also explained by fewer conflicts among aircraft leading to a fewer number of edges with infeasible time windows. The lower number of conflicts is evident from the sum of columns I 1 and I 0 which indicate the number of potential conflicts on the unimpeded route. The conflicted aircraft range from 4% to 14% in instances marked original. As a result, speed profiles in the multi-graph with u = 1 can comfortably satisfy time window constraints in most cases.
Secondly, we analyse the instances marked ins1 and ins2. These instances denote artificial scenarios created as detailed in "Methods" section which introduce conflicts among aircraft on purpose. The higher proportion of conflicts in the instances resulted in larger relative savings up to around 3% with higher values of u, particularly for the HKG airport. The savings with u = 5 are higher than for u = 3 . In most cases, ins2 instances obtained less savings compared to ins1. Instances in ins2 are similar to ins1 but with less conflicts caused by sequential routing of aircraft. The sequential routing can lead to both better or worse routes for the subsequent aircraft. The difference in results for ins1 and ins2 reflects the effect of sequential routing. As can be seen, even without Lastly, artificial instances with increased traffic are investigated. These instances are marked 0-100 denoting the % increase in traffic over the original levels at the rush hour. When the traffic is high, several aircraft routes are likely to be in conflict with each other. The results show mostly positive values of savings which are increasing with higher traffic. As with previous experiments, the savings with u = 5 are higher than for u = 3 . The trend in improvement with higher levels of traffic is clearly documented in the increasing number of routes denoted as 'improved' . As the traffic increased, the number of conflicts and potential routes which can be improved with more speed profiles are increased too. For DOH and PEK airports, the % savings are both positive and negative. The DOH instances have less conflicts (as indicated in the values of I 1 and I 0 ) compared to other airports. Also, the airport layout is simple and the routing algorithm has less options to take a detour. Instead, the routing algorithm with u = 1 delays the start time of an aircraft where the extra waiting time does not constitute additional fuel costs as it happens at the start/end of the route. For departures this is achieved at the gate with engines turned off. For arrivals, it is assumed that postponing can be achieved before landing via air traffic control procedures and therefore not affecting fuel costs at ground. With a higher u, the route can be found without postponing but with higher costs than u = 1 and extension. For PEK airport, the magnitude of savings and improved routes are both lower than for HKG. This can be explained by a more complicated layout of PEK than HKG. As a result, the opportunity of taking a detour is more frequent and less penalising (shorter detours possible) than at HKG. This causes less savings for the routes without detours found with higher u.
The indicated number of routes which can be improved in column I 1 is higher than the number of routes which are actually improved. The ratio I eff which is the number of 'improved' routes to I 1 when u = 5 ranges from around 13% to 87%. The ratio of I 1 to the total number of aircraft in the instance ranges 4-50%. This means, that if the proposed indicator is used only 4-50% of original aircraft have an opportunity for improvement with higher values of u. Table 4 used the multi-graph with nondominated speed profiles. This makes sense if there are no time window constraints. Figure 5a shows an objective space for a segment. If we assume that the Pareto front for speed profiles is continuous then any solution can be generated on the indicated curve. The assumption comes from the fact that the variables for the speed profile generation are continuous 14 . The curve has usually a parabolic shape and the left part constitutes the Pareto front. The right part has an upward trend, as speed profiles with very long taxi time after some threshold consume more fuel than faster speed profiles. For any feasible dominated solution (solution a), there is a corresponding nondominated solution (solution b) with the same taxi time and better fuel consumption. Therefore, to reach any solution on the feasible part of the Pareto front, only nondominated solutions located on this front have to be included in the multi-graph. However, when the infeasible region covers a larger portion of the Pareto front and extends beyond the last nondominated solution (solution b) as in Fig. 5b, the previously dominated speed profile (solution c) becomes nondominated. Therefore, to reach this new part of the Pareto front, the previously dominated speed profile needs to be included in the multi-graph.

Dominated speed profiles. The experiments presented in
The presence of time window constraints affects the search in AMOA*. If some edge does not have any time window available for the current speed profile, then holding can be applied (at additional cost) or another speed profile can be selected if available. In the following, we conducted experiments on a multi-graph with u = 1 , u = 3 and u = 5 with holding where the aircraft is held at the end of the previous segment with engines running at idle until the time window becomes feasible. Also, in another set of experiments, instead of holding, the routing algorithm checks the database of dominated speed profiles and selects one with the least cost which is feasible. Both approaches are illustrated in Fig. 6. It should be noted that in both cases we are effectively adding a dominated speed profile into the multi-graph. Table 5 details cost savings obtained with applied holding and dominated solutions. The column headings (e.g. cost(u = 3 )) denote the baseline case without holding against which the results are compared with. Columns

Discussion
In this study, we analysed cases and conditions when the multi-graph structure directly influences the quality of the obtained solutions found by a routing algorithm. The airport ground movement problem was employed as a case study. The experiments were carried out with the different number of the parallel edges and costs. The higher number of parallel edges resulted in higher number of Pareto optimal solutions found and time complexity of the search. Different approaches for speed profile selection could control which parts of the Pareto front are covered by the solutions.
For the original instances, the increased u resulted in negative or only small improvements < 1% . However, artificial instances with higher number of conflicts showed higher savings up to 4.12% in some cases. With higher traffic, the savings increased too. This result suggests that when traffic levels are low the multi-graph with u = 1 can find good quality solutions in most cases. As search with u = 1 is fast compared to a case with higher u, using u = 1 is preferential. Therefore, the routing algorithm can be adopted for a multi-graph with u = 1 . In this case, the multi-graph can be constructed using preferences w p determined by the decision maker, or w p can be iteratively changed to cover different parts of the Pareto front. For the higher traffic levels, u = 3 can improve the results and u = 5 even more.
An indicator proposed to estimate the number of routes which can be potentially improved by using alternative speed profiles could be successfully used to indicate in which case a higher value of u is beneficial. Therefore, the routing algorithm should use a higher values of u only for those instances (i.e. aircraft), while for others, u = 1 is sufficient. The ratio of I 1 to the total number of aircraft in the instance ranged 4-50%. If u > 1 is used only for those aircraft indicated in I 1 , 50-96% aircraft can use u = 1 and save computational time.
Also, the results highlighted the differences in savings between different airports, pointing to the importance of the airport layout. On the other hand, as savings with higher u demonstrate, using more complex modelling techniques can bring benefits in high traffic scenarios and potentially offset costly investments in new taxiways.
The experiments with dominated speed profiles showed the importance of considering including dominated parallel edges in a multi-graph in the presence of time window constraints. With dominated speed profiles included, the routing algorithm could find routes up to 7.85% better compared to the multi-graph with the same u without dominated speed profiles.
In future work, a similar case study could shed light on the creation of the multi-graph for other problems based on multi-graphs such as the vehicle routing problem, hazardous material transportation, the multimodal shortest path problem and trajectory based traffic management. Also, the routing algorithm could actively use the   Cost(u = 1) Improved Cost(u = 3) Cost(u = 1) Improved Cost(u = 5) Cost(u = 1) Improved www.nature.com/scientificreports/ indicator or some similar measure to create multi-graphs with a structure which is fast for search but guarantees high quality solutions. The results in this study also highlighted the issue of sequential routing when the route of one aircraft can affect the subsequent aircraft. One way of addressing this problem is to search for a better sequence of aircraft. However, even more promissing global approach would be considering multiple aircraft simultaneously such that the multi-graph structure could be utilised to avoid conflicts not only with subsequent aircraft but all aircraft considered.

Methods
Multi-graph. The airport layout is represented as a directed graph G = (V , E) . Nodes n ∈ V represent gates, stands, taxiway intersections, intermediate points and runway exits. Edges e ∈ E represent taxiways between two nodes. A sequence of edges of a similar type between two nodes n, m is defined as a segment s n,m = (e 1 , e 2 , . . . , e h ) . The segments are of two types, straight and turning. If an edge and its predecessor edge (in the direction of a taxiing aircraft) have an angle ≥ 30 degrees 15,16 , then it will belong to a turning segment. Otherwise, it is part of a straight segment. Consecutive edges of a similar type (straight, turning) are grouped together. For a segment between two nodes n, m, (c n,m,l,1 , c n,m,l,2 , . . . , c n,m,l,q ) ∈ C n,m is a cost vector with q objectives which corresponds to the lth speed profile for that segment. Speed profiles for a single segment are continuous functions of time 14 .
In this study, a piece-wise linear function with four phases 14 including acceleration, constant speed, deceleration and rapid deceleration is adopted. The duration of each phase and the associated thrust levels determine the taxi time and fuel consumption of the speed profile. Evenly distributed speed profiles according to the two objectives from Pareto solutions for all segments were adopted in our previous study 17 . The speed profiles include different weight categories of aircraft for up to 20 speed profiles per segment. The multi-graph is constructed using these speed profiles.
As multi-graphs with a large number of speed profiles per segment can be computationally prohibitive, two different approaches are proposed to reduce the multi-graph: (1) From the Pareto front of speed profiles for a segment, u evenly distributed solutions can be considered. In this study, rows in C n,m were ordered according to the first objective and u solutions were selected with even distance from each other according to the first objective. (2) If preferences for the search are known beforehand, the first u speed profiles ranked according to that preference can be considered. In this study, we take the scalar product of obj 1 and obj 2 in each row in C n,m and unit costs w P as the preferences. These aggregated costs then give the ranking of the rows in C n,m and lower costs are preferred. As an example, consider a Pareto front of speed profiles with two objectives: {(1,6),(2,4),(4,3),(5,2),(7,1)} where we want to select u = 3 solutions. In the case of evenly distributed solutions, we select the extreme solutions (1,6) and (7,1) and the middle solution (4,3) which has the same distance from the extreme solutions w.r.t. the first objective. This way, the Pareto front is evenly covered with our selection. In the case of known preferences, e.g. w p = (M, 0) where M is a large number, we select (1,6) with the lowest cost and two solutions (2,4) and (4,3) with the 2nd and 3rd lowest cost, respectively. Airport layout. The taxiway layout was processed 18 from OpenStreetMap (www.openstreetmap.org). All edges were set as bi-directional. The flights for DOH and HKG instances were captured with the specialised tools 18 from freely-available data on FlightRadar24.com. The data specifies landing/pushback times, gates/runway exits and weight category for each flight.
The instances marked ins1 and ins2 are artificial instances. From the original instances, pairs of aircraft with overlapping unimpeded routes were selected and their start time changed such that they arrived to where their routes intersected at the same time . When the traffic is high, several aircraft routes can be in conflict with each other. In order to eliminate such interference, pairs of aircraft from ins2 are separated from other pairs by a large time interval.
The instances marked 0-100 were generated as follows. From the original instances, 1 hr of traffic during rush hour was selected, marked 0 (e.g. doh0). Then, artificial instances were generated with traffic levels increased by 25, 50, 75 and 100% by randomly adding additional aircraft. It should be noted that some instances would be unrealistic due to the fact that the runway is usually the main bottleneck at airports. With a theoretical maximum capacity of 60 aircraft per hour for a single runway, some instances for HKG (2 runways in total) and DOH (1 runway) could exceed this capacity. Improvement indicator. An indicator is proposed in this paper to indicate the number of routes that could be improved by using alternative speed profiles. The indicator is outlined in Algorithm 1. For an aircraft, the routing algorithm is run in Line 1 to find and select a route (using the same w P = (0.469, 0.71) as above) without considering any time windows (the unimpeded route). This step can be carried out offline before the search. Then, during the search the routing algorithm is run again in Line 2 considering time windows. The resulting route is compared with the unimpeded route in Line 3. If the two routes are identical, then the route is already the best route and cannot be improved by using different speed profiles. Otherwise, the difference is caused by a conflict with another aircraft and the corresponding time window causing the conflict is found. Firstly, the first differing edge (in the direction from the start to destination node) is identified in Line 7. Then, the fastest speed profiles are applied in Line 8 to calculate the fastest time t 1 at which the aircraft can arrive at the identified edge. Furthermore, t 2 is computed as the time of applying the slowest speed possible (5.14 m/s, 10 knots) in Line 9. The interval between t 1 and t 2 corresponds to times at which the aircraft can arrive at the edge using different speed profiles on the previous segments. If this interval in Line 13 is wide enough to traverse the edge using minimum www.nature.com/scientificreports/ traversing time (the time using the fastest speed profile) and not in conflict with the time window, then the route is included in I 1 in Line 14. The route counted in I 1 can be potentially improved by avoiding the conflicting time window using different speed profiles. Otherwise, the route is included in I 0 in Line 16 because no matter which speed profile is used, the arrival and departure times from the edge will conflict with the time window. I 0 is therefore the number of routes which cannot be improved even by using alternative speed profiles. It should be noted, that I 1 only indicates a potential improvement and serves as a metric corresponding to the upper bound on how much the results can be improved using u > 1 speed profiles. The alternative speed profile can be more costly than the detour and also any time window conflicts on subsequent edges are not tested. I eff is the ratio of actually improved routes to I 1 using u > 1 speed profiles. The actually improved routes used alternative speed profiles to avoid the conflicting time window and the subsequent detour which resulted in better objectives.
It should be noted, that the proposed indicator does not compromise the search. If there are no conflicting time windows, the selected route from the Pareto front found by the routing algorithm with u = 1 and u > 1 are the same due to the minimum cost being the same. This condition is tested in Line 3. Otherwise, the indicator algorithm tests if it is feasible to prevent the conflicting time window. In the case of I 1 , the conflict can be prevented with u > 1 and full search should be conducted. In the case of I 0 , the conflict is inevitable and therefore search with u > 1 cannot find a better route.

Data availability
The airport layout data sets and anonymised aircraft traffic instances are available here: https:// github. com/ mweis zer/ amoa_ repor ts/ blob/ main/ data. zip.